function val = FOC(x)

global parms lnVcoef DlnVcoef stateV
% parms = [delta A Sbar alpha0 alpha1 alpha2 alpha3 Gamma1 alpha psi beta gamma Q0 popgr Abar];
% % Note that parms is a 1x15 vector

delta = parms(1);
A = parms(2);
Gamma1 = parms(8);
alpha = parms(9);
psi = parms(10);
beta = parms(11);
gamma = parms(12);

k = stateV(1);
H = stateV(2);

inv = x(1);
j = x(2);

B = A*k;
psc = 1-psi;
c = (1-(Gamma1+H).^(-alpha))*B-inv-j;

knxt = inv+(1-delta)*k;
Hnxt = H+B^psi*j^psc;

lnk = log(knxt);
lnH = log(Hnxt);
var_poly = [1,lnk,lnH,lnk.^2,lnH.^2,lnk.*lnH];
lnV = var_poly*lnVcoef';
V = -exp(lnV);

Dvar_poly = [1,lnk,lnH];
DlnV = Dvar_poly*DlnVcoef;
DlnVk = DlnV(:,1);
DlnVH = DlnV(:,2);

val = zeros(2,1);

val(1) = -c^(-gamma)+beta*V*DlnVk/knxt;
val(2) = -c^(-gamma)+beta*(1-psi)*V*DlnVH*B^psi*j^(-psi)/Hnxt;
% val(2) = DlnVH*(1-psi)*B^psi*j^(-psi)/Hnxt-DlnVk/knxt;

return